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ABSTRACT 

Weak interstellar scintillations of pulsar B0809+74 were observed at two epochs 
using a 30m EISCAT antenna at 933 MHz. These have been used to constrain 
the spectrum, the distribution and the transverse velocity of the scattering plasma 
with respect to the local standard of rest (LSR). The Kolmogorov power law is a 
satisfactory model for the electron density spectrum at scales between 2 x 10 7 m and 
10 9 m. We compare the observations with model calculations from weak scintillation 
theory and the known transverse velocities of the pulsar and the Earth. The simplest 
model is that the scattering is uniformly distributed along the 310 pc line of sight 
(/ = 140°, b = 32°) and is stationary in the LSR. With the scattering measure as the 
only free parameter, this model fits the data within the errors and a range of about 
±10 km s _1 in velocity is also allowed. The integrated level of turbulence is low, 
being comparable to that found toward PSR B0950+08, and suggests a region of low 
local turbulence over as much as 90° in longitude including the galactic anti-center. 
If, on the other hand, the scattering occurs in a compact region, the observed time 
scales require a specific velocity-distance relation. In particular, enhanced scattering 
in a shell at the edge of the local bubble, proposed by Bhat et al. (1998), near 72 pc 
toward the pulsar, must be moving at about ~ 17 km s -1 ; however, the low scattering 
measure argues against a shell of enhanced scattering in this direction. The analysis 
also excludes scattering in the termination shock of the solar wind or in a nebula 
associated with the pulsar. 

Subject headings: scattering — plasmas — ISM: kinematics and dynamics — turbu- 
lence — radio continuum: ISM — pulsars: individual (B0809+74) 
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1. Introduction 

The use of pulsar scintillation to probe the in- 
terstellar medium has been complicated by the fact 
that most observations below about 1 GHz are in 
the strong scintillation regime. In this regime there 
are both refractive and diffractive components in the 
intensity scintillation pattern, and neither spatial 
scale is known a priori. Furthermore the contribu- 
tions from various locations along the path of prop- 
agation do not add in a simple linear fashion (see 
Rickett, 1990 and Narayan, 1992 for reviews). Weak 
scintillation is much simpler to interpret since con- 
tributions from various locations add linearly and 
the spatial scale of each contribution depends only 
weakly on unknown parameters such as the spectral 
exponent and the inner scale. The temporal scale of 
each contribution depends primarily on its distance 
and its velocity In principle, it is feasible to "invert" 
the observations to estimate the spatial spectrum of 
the electron density and, by exploiting the variable 
velocity of the Earth, to investigate distribution of 
the scattering plasma and its velocity. 

Observations can be made in weak scintillation 
by using higher frequencies and/or nearby pulsars. 
However it is difficult to obtain a stable estimate 
of the covariance function because the time scale of 
the scintillation is of the order of an hour. Thus a 
given observation will typically contain only a few 
time scales. Earlier observations by Backer (1975) 
at 3 GHz and Malofeev et al (1996) at 3, 5 and 
8 GHz have confirmed that weak scattering occurs 
more or less as expected, but have not been long 
enough to provide the statistical accuracy necessary 
for detailed model fitting. 

In this paper we report two longer measurements 
of a circumpolar pulsar in weak scintillation. The 
temporal statistics of these data sets are adequate to 
estimate the spectrum over a decade in wavenumber 
and put constraints on the distribution and velocity 
of scattering plasma. We have developed a proce- 
dure for the analysis that may be more generally 
useful. 

2. Observations and Data Analysis 

The pulsar B0809+74 was selected because ear- 
lier observations at lower frequencies suggested that 
it should be in weak scintillation at 933 MHz, which 
is the only frequency available on the 30m EISCAT 
antenna at Sodankyla, Finland, and because it is 
circumpolar at that latitude (68N). The receiver 



passband is 8 MHz centered on 933 MHz and one 
linear polarization was used. The use of one linear 
polarization introduces a potential bias in the analy- 
sis since the pulsar is partially polarized. As a result 
the apparent flux varies due to the rotation of the 
antenna and the Faraday rotation in the ionosphere. 
This was thought to be negligible at 933 MHz, but 
was later found to be important. The effect was 
eliminated by selecting an unpolarized portion of 
the pulse profile as discussed below. 

Impulsive interference is common at Sodankyla 
and the square law detector was designed to reduce 
the contribution from very short noise pulses. The 
pulsar intensity could be adequately sampled at 100 
samples per sec (sps), however it was actually sam- 
pled at 100 Ksps using a very broad band detector. 
Pulses exceeding about 10 standard deviations were 
clipped at the sampler. The 100 Ksps series was 
then smoothed and decimated by a factor of 1000. 
This process reduced the contribution of impulsive 
interference to a low level without distorting the pul- 
sar signal, and also reduced the rms quantization 
noise by a factor of 30. We observed continuously 
for 75 hours starting at 20 hr UT on 1996 April 
8, although the observations were stopped due to 
system errors and restarted several times. We also 
observed a second epoch from 23 hr UT on 1998 
September 30 for 85 hrs. 

The data were summed (off-line) according to the 
apparent pulsar period (Doppler-shifted to the ob- 
servatory). Each period (ssl.29 s) was tested for 
interference using the ratio of the rms to the mean 
off-pulse power. The expected value of this ratio is 
about 0.003 and its standard deviation for a single 
pulse period is about 0.0003. Periods for which this 
ratio exceeded 0.005, about 7 standard deviations, 
were rejected. Successive integrations of the pulse 
profile were computed every 2 minutes. These were 
examined by hand and obvious interference was re- 
moved. The automatic interference test was very 
effective so less than 1% of the 2 minute averages 
were edited manually. Pulse profiles were also com- 
puted every 10 minutes. Both 2 and 10 minute av- 
erages were used in subsequent work. The average 
pulse profile was used to estimate the energy in an 
on-pulse window and also in an equal off-pulse win- 
dow. First an estimate of the background power 
was subtracted, then the energy in each of the two 
windows was computed. 

Initially we used an on-pulse window which in- 
cluded the entire pulse, however a referee, (Michael 
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Kramer), pointed out that B0809+74 has 30% linear 
polarization (Gould and Lyne, 1998). Examination 
of our pulse profiles revealed a 24 hr periodicity in 
the leading edge, which could be attributed to the 
linear polarization of the pulsar, the rotation of the 
telescope, and the diurnal change in the Faraday 
rotation in the ionosphere. This contributed about 
10% of the variance in the pulse flux. Although 
this is not large it was greater than our estimation 
error and significantly altered our conclusions. The 
Gould and Lyne observations show strong linear po- 
larization in the leading half of the pulse with no ro- 
tation of position angle, but the trailing half of the 
pulse is essentially unpolarized. Thus we changed 
our pulse energy computation to use a 50 msec on- 
pulse window starting at the pulse peak and includ- 
ing all the trailing half. With this window the polar- 
ization modulation was not detectable. The result- 
ing time series of on-pulse power and off-pulse noise 
are displayed in Figure |l| with 10 minute resolution. 
The flux density was not calibrated absolutely, it 
has been normalized to unity averaged over the en- 
tire observing period. 

In the absence of explicit flux density calibration, 
we used the variation of the system noise level to ob- 
tain a limit on the gain variation. The off pulse sys- 
tem noise varied almost regularly with a period of 24 
hrs and an amplitude of about 6% peak to peak. We 
assume that most of this variation is due to changes 
in the atmospheric and ground noise contributions, 
as the zenith angle varied from 6° to 38° and the 
azimuth remained in the northern quadrant during 
the 24 hr period. We therefore conclude that the 
gain variation is periodic over 24 hr with an ampli- 
tude considerably less than 6% and thus is negligible 
compared to the large variations over a few hours 
observed in the pulse intensity. 

The second order temporal statistics, i.e. the 
power spectrum and the autocovariance function, 
can be theoretically modelled and used to estimate 
the distribution of scattering material, its velocity 
and some characteristics of the turbulent spectrum. 
Although the two statistics are simply related by a 
Fourier transform, their estimation errors are quite 
different and both are useful for different purposes. 

To estimate the time scale we computed the 
structure function rather than the autocovariance 
function, since the structure function is a more re- 
liable estimator when the number of independent 
samples of the scintillations is not very large, as is 
the case here. For this purpose we used the 10 min 
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Fig. 1. — Observations of the total pulse flux of 
B0809+74 at 933 MHz, normalized to the mean flux 
over the observing period. The samples plotted are 
10 minute averages. The lower traces are computed 
from an off-pulse window to provide a noise esti- 
mate, (a) Data from April 1996. (b) data from 
October 1998. 

data because the noise correction for these data is 
negligible and the resolution is more than adequate. 
The structure function was computed in the usual 
way from the series of pulsar energies Ij (normalized 
by their estimated mean) sampled at times tj. 

D^^^^-i.r/^i (i) 

Here Ej J represents the summation over all data 
pairs with time difference ti — tj lying within a win- 
dow of width Sp centered on 77 = ISp, where I is an 
integer and 5p is the resolution of the data. The 
maximum lag computed was three quarters of the 
total duration of the observations. 

Errors in the structure function result from ad- 
ditive receiver (radiometer) noise, intrinsic pulse to 
pulse variation, and from estimation errors due to 
the finite observing span. The receiver noise and in- 
trinsic variation are effectively white and indepen- 
dent of the scintillation. They represent additive 
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constants in the structure function which can be 
subtracted. The receiver noise contribution can be 
estimated very accurately from the off-pulse noise. 
To measure the intrinsic variation we used the 2 
min data. The first point in the structure function 
of the 2 min data is dominated by intrinsic varia- 
tion and receiver noise. This was used to estimate 
the contribution of intrinsic variation to the struc- 
ture function of the 10 min data with very good 
accuracy. The corrected structure function of the 
10 min data is shown in Figure 0. 
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Fig. 2. — Structure functions of pulse flux time se- 
ries shown in Figure [l] corrected for receiver noise 
and intrinsic variation. The data are marked by 
±a error bars: (a) April 1996; (b) October 1998. 
The solid lines are the best fit theoretical models 
for a thin screen. The dashed lines are best fit mod- 
els for a uniform extended scattering medium (as 
described in section 3.2). The dots represent the 
structure function of the receiver noise. 

The estimation errors in the structure function 
due to the finite data length are substantial; they 
are derived in Appendix B. The ±er estimation er- 
rors, given by equation ( B12 ), are indicated by the 
error bars in Figure ||. As discussed in the Appendix 



the fractional error due to estimation error decreases 
toward small time lags but is partially compensated 
by the addition of the errors due to noise. The esti- 
mation errors are correlated over about an octave in 
time lag below saturation. At larger lags they are 
correlated over a characteristic time scale defined 
as the 50% width of the autocovariance. For model 
fitting purposes, the lags were decimated to reflect 
this distribution of independent estimates. 

The power spectra were computed from the 2 min 
data using a direct Fourier transform procedure. 
The 1996 data was divided at the data gap and the 
two halves were transformed separately and their 
spectra added. The 1998 data was transformed in 
one block. The spectra were computed with a raw 
resolution of 7.3 xlO -3 cycles per hour (cph) and 
normalized so that the integral under the spectrum 
is the variance. They were boxcar smoothed and 
decimated by a factor which increased from 3 for 
the lowest 4 points to 125 for the highest frequen- 
cies. This reduces the estimation error and broadens 
the resolution. The white noise contribution from 
receiver noise and instrinsic variation was estimated 
from the frequencies above 5 cph, where the scintil- 
lation power is negligible, and subtracted from the 
spectral estimate. The final power spectra are plot- 
ted in Figure 3. The vertical bars indicate the ±cr 
estimate error and the horizontal bars indicate the 
frequency resolution. Power spectra are convenient 
for model fitting, particularly at the high frequen- 
cies, because the estimation errors are independent 
and their variance is known. Thus it is easy to ob- 
tain reliable confidence intervals. 

Two basic parameters which describe the flux 
time series are the normalized standard deviation 
or scintillation index m w , and the time scale r w at 
which the correlation function falls to 50%. They 
were derived from the structure function using the 
fact that Di(t) saturates at 2m^, for large time lags 
and estimating the time scale as the lag where Dj(t) 
crosses half the saturation level (i.e. Dj{t w ) = m^). 
By this method we found m w = 0.77 ± 0.2 and 
t w = 1.93 ± 0.2 hr for the 1996 observations and 
m w = 0.84 ± 0.15, and r w = 2.19 ± 0.2 hr for the 
1998 observations. 

3. Model Fitting and Interpretation 

Even though the observed scintillation indices 
(« 0.8) are not much less than unity, we initially in- 
terpreted the results using weak scintillation theory 
The advantage of weak scintillation theory is that 



4 




0.01 0.1 1 10 



(a) Frequency (cycles/hr) in 1996 




0.01 0.1 1 10 
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Fig. 3. — Power spectra of the pulse flux time se- 
ries shown in Figure [l] corrected for receiver noise 
and intrinsic variation, (a) 1996 and (b) 1998. The 
data are marked by ±<r error bars and the resolu- 
tion is indicated by the horizontal bar. The off- 
pulse noise spectrum is marked by a dotted hori- 
zonal line. The sum of the off-pulse noise and the 
intrinsic noise is marked by a dashed horizontal line. 
The frequency range from which the intrinsic contri- 
bution was estimated is marked with a heavy hor- 
izontal line. The best fit theoretical models for a 
thin screen are shown by solid lines. The best fit 
theoretical models for a uniform medium which is 
stationary in the LSR are shown by dashed lines. 
The best fit models in 1996 have zero inner scale, 
but in 1998 they have s, = 2x 10 7 m. A best-fit 
theoretical model for a uniform medium with zero 
inner scale is also drawn over the 1998 data as a 
dotted line for comparison. 

the observed structure functions and power spec- 
tra are linear integrals contributed from all points 
on the line of sight as described by the Born ap- 
proximation discussed in Appendix A. Under weak 
scintillation we can do quantitative modelling of the 
density spectrum and the velocity distribution of 
the scattering medium. However when the scintilla- 



tion index is approaches unity, there can be signif- 
icant departures from the linear Born model. Sim- 
ulations by Frehlich (1999, private communication) 
corresponding to our observed m w « 0.8 show that 
the primary effect of incipient strong scattering is to 
broaden the spectrum by a factor of 1.20±0.02 and 
to reduce the time scale t w by the same factor. The 
scintillation index m w is also reduced by 4% from 
the weak scattering approximation. The statistical 
error in m w is about 20% so the 4% bias due to 
incipient strong scintillation is negligible. However 
the statistical error in t w is only 10%, so we have 
divided all our Born model calculations by the fac- 
tor of 1.20 to ensure that the bias is less than the 
statistical error. Strong scattering also smoothes 
the "Fresnel ripples" in the spectrum and rounds 
the "knee" of the structure function, but this is a 
second order effect. 

We assume that the density spectrum is an isotropic 
Kolmogorov power law, and also consider the effect 
of an inner scale. We need only consider one basic 
geometry: the pulsar is a point source at a distance 
L from the Earth and a thin scattering screen is 
located z p from the pulsar and z a from the earth, 
where z p + z a = L. To compute the temporal statis- 
tics we also need the velocity of the Earth, the pul- 
sar, and the scattering screen. We can then com- 
pute the effect of a distribution of scattering ma- 
terial along the line of sight by integrating over all 
such elementary screens. 

We have referred all velocities to the "local stan- 
dard of rest" (LSR), since that appears to be the 
logical reference frame for the velocity of the inter- 
stellar plasma between the Earth and the pulsar. In 
this frame the solar velocity is 20 km s _1 directed 
towards right ascension of 18 hrs and declination of 
+30°. The pulsar velocity (V a = 31 ± 10 km s -1 ; 
Vs = —55 ± 9 km s" 1 in this frame) is determined 
from its proper motion (Lyne, Anderson, and Salter, 
1982) and its distance (310 pc) from the model of 
Taylor and Cordes, 1993 [TC93]. 

In the following we will compare the observa- 
tions with two models: (1) a single thin screen and 
(2) a uniform distribution. We find satisfactory 
agreement for both models, and then explore what 
range of LSR velocities is allowed for the scattering 
medium. 

3.1. Thin Screen Model 

The intensity structure function for weak scat- 
tering by a thin screen with a simple Kolmogorov 
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spectrum (equation A12), has the universal form 
(A13). The time scale r w = 0.98 rf/V c g, where 
77 = yJz z p /Lk is the "Fresnel scale" , i.e. the spa- 
tial scale at the screen that contributes most of the 
scintillation, and V e g is the velocity of the scatter- 
ing plasma across the line of sight given by equation 
( |A9| ). (It should be noted that the spatial scale at 
the observing plane is expanded by spherical diver- 
gence to r/L/zp). 

The best fit thin screen models are drawn over 
the observed structure functions as solid lines on 
Figure |[ The structure function at the smallest lag 
was excluded from the fit, in order to avoid the pos- 
sible influence of an inner scale and because the cor- 
rection for intrinsic pulse variation has uncertainties 
not reflected in the error bars at the smallest lag. As 
discussed earlier the lag values were chosen to sam- 
ple approximately independent estimates and two 
parameters were fitted. This gave approximately 
28 and 34 degrees of freedom, respectively in 1996 
and 1998. The minimum \ 2 were 34 and 30, re- 
spectively, which indicates satisfactory fits to the 
data; the differences reflect the fact that the errors 
in the saturated portion of the structure functions 
are noticeably different in the two years. We use 
these fitted models to make refined estimates of the 
time scale and its errors: t w = 1.99 ± 0.18 hr and 
r w = 2.05 ± 0.16 hr for the 1996 and 1998 obser- 
vations, respectively. These time scale estimates 
average over a range of time lags and so reduce 
the effect of the particular errors near where the 
structure function crosses half its saturation level. 
Hence we also use these same time scale estimates to 
characterize the data in section 3.2 on the extended 
medium model. 

Weak scintillation theory predicts the relative 
contribution to the variance m^, and the Fresnel 
scale Tf for a screen at position < z Q < L, which 
are plotted in panels (a) and (c) of Figure ||[ The 
effective velocity V e g for the mean time of the 1996 
observations is plotted as a solid line and that for 
the 1998 observations is plotted as a dashed line in 
panel (b) under the assumption that the scattering 
plasma is stationary with respect to the LSR. The 
predicted time scales = 0.82 rf/V e s (corrected 
for m w =0.8) are plotted in the lower panel Figure 
||(d) as a thick solid line for 1996 and a thick dashed 
line for 1998. The uncertainty limits in these time 
scales caused by errors in the pulsar proper motion 
are indicated by light solid and dashed lines at ±<r. 
Also shown are the confidence intervals for the mea- 
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Fig. 4. — Dependence of weak scintillation param- 
eters on screen distance: (a) contribution to vari- 
ance; (b) effective velocity V c g with respect to the 
LSR (solid line 1996, dashed 1998); (c) Fresnel scale; 
and (d) time scale (solid line 1996, dashed 1998). 
The time scales have been reduced by a factor 1.20 
from the Born theory to account for incipient strong 
scattering. The thin lines show the envelope of time 
scales corresponding to errors in the pulsar proper 
motion. The confidence intervals for the observed 
time scale are shown as horizonal line pairs near 2 
hrs (solid line 1996, dashed 1998). The bars indi- 
cate allowed distances to the screen at each epoch 
assuming the nominal pulsar velocity. The rectan- 
gle indicates the distance to the shell of enhanced 
scattering proposed by Bhat et al. (1998). 

sured time scales as horizontal pairs of lines. 

Neglecting the error in the proper motion for the 
moment, one can see that there are two locations 
at which a thin screen would provide the time scale 
measured in 1996, at 51 ± 6 pc and 213 ± 9 pc. 
After analyzing the 1996 observations we realized 
that this ambiguity could be resolved with addi- 
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tional observations at a different epoch (we thank 
Don Backer for pointing this out). Then the Earth's 
velocity would be different and so curves for V e s(z ) 
and t w (z ) would also be different. This is shown 
in Figure || where we plot V e e(z ) and r w (z G ) for six 
observing epoch's spaced every two months on the 
8th of each numbered month. One can see that the 
velocity near the pulsar does not change much but 
the velocity near the Earth changes dramatically. If 
the scattering screen were located near 213 pc, then 
the time scale would remain near 2 hrs. However if 
the scattering screen were located near 51 pc, the 
time scale would be much larger if observations were 
made between July and October. Accordingly we 
chose to reobserve the source in October 1998. One 
can see from Figure ||(d) that the 1998 observations 
are consistent with distances of 1±0.1 pc or 178±12 
pc. The ambiguity is thus resolved in favor of the 
larger distance. However the match is not perfect as 
the ±<7 error bars do not overlap. The probability 
of this occuring by chance is less than 2% so we are 
led to consider the possibility that the medium may 
not be stationary in the LSR, which we discuss later 
in Section 4. We note that the uncertainty in screen 
location due to errors in the pulsar proper motion 
are larger than those due to errors in the measured 
time scale. Thus the two possible screen distances 
(213 ± 9 pc and 178±12 pc) are in fact jointly un- 
certain by about ±28 pc. However, since the pulsar 
velocity must be the same at each epoch there is 
no single distance that matches both observations. 
It should also be noted that, in any case, it is very 
unlikely that the 1998 observations could have come 
from the closer location, because that would require 
a very thin, very intense scattering screen. 

We also fit the simple Kolmogorov model with 
no inner scale to the estimated power spectra. The 
1998 data were best fit with a nonzero inner scale 
(as described in section 3.3) but this did not af- 
fect the estimate of the time scale. The results, 
r w = 1.84 ± 0.22 hr for 1996 and r w = 1.84 ± 0.15 
hr for 1998, are not quite consistent with the time 
scales derived from the structure functions. The 
confidence limits for the two methods just touch, 
suggesting that there may be a real bias. This could 
occur because the power spectra and the structure 
functions emphasize different aspects of the data. 
The structure function fit is dominated by scales 
near t w whereas the power spectrum fit is domi- 
nated by the higher frequencies. The effect is that 
estimates from the power spectrum would move the 
location of the screen about 15 pc closer to the pul- 
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Fig. 5. — Dependence of weak scintillation param- 
eters on screen distance plotted for six epochs as 
indicated by month/day: (a) V e g with respect to 
the LSR; (b) Scintillation time scale (with 20% re- 
duction to correct for incipient strong scattering). 

sar. We use the time scales from the structure func- 
tion in further analysis. 

The plots of Figures | and | are based on models 
in which the screen is at rest in the LSR. In the 
Discussion (section 4) we consider models in which 
the screen is moving with respect to the LSR, and 
compare these with other astronomical information 
on motion in the local ISM. We now consider models 
with a uniform distribution of scattering. 

3.2. Extended Scattering Model 

When the scattering medium extends from the 
source to the observer, the theoretical structure 
function is the sum of many screen contributions 
each with an effective velocity V e s(z ) as in Ap- 
pendix A. One can visualize this by looking at Fig- 
ure |J we add screen structure functions weighted 
by the curve in the top panel, with each structure 
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function scaled in time by the heavy curves in the 
bottom panel. We computed theoretical structure 
function for the 1996 and 1998 observing times as- 
suming that the scattering plasma is uniformly dis- 
tributed along the line of sight and is stationary in 
the LSR. As noted above we reduced the time lags 
by a factor of 1.20 in each model, to correct for 
the effect of incipient strong scattering. Then with 
the velocities known, the total variances m%, at each 
epoch are the only free parameters in the models. 
The corresponding least squares fits are plotted as 
dashed lines in Figures ||a and b. These models 
have a more rounded transition to saturation than 
for the screen geometry, and they fit the observed 
time scales reasonably well at half the saturation 
level. The \ 2 values for these fits are 39 and 30 
for 1996 and 1998 with 28 and 34 degrees of free- 
dom, respectively. This shows that the extended 
medium model fits the data as well as the best fit 
thin screen and is intuitively more appealing as it 
has only one free parameter. We note that incipient 
strong scattering causes a slight change in shape of 
the structure function as well as in time scale; hence 
we cannot use the differences in \ 2 to discriminate 
between the models. 

The LSR is based on the motion of nearby stars 
and it may not be the best reference frame for the 
scattering plasma. Though differential Galactic ro- 
tation can be neglected over the 310 pc to the pulsar, 
we have to consider what is known about motion of 
the nearby interstellar gas. Current views of the "lo- 
cal bubble" are described in the proceedings of IAU 
colloquium 166 (Breitschwerdt et al 1998). This ion- 
ized region within 100-200pc of the Sun must largely 
determine the radio scattering for our pulsar obser- 
vations. However, the basic idea of a hot ionized 
low density cavity has been refined by various recent 
measurements. Lallement (1998) discusses measure- 
ments of local interstellar "cloudlets" within a few 
pc of the Sun moving at 5-10 km s^ 1 with respect 
to the LSR. Genova et al. (1998) describe the kine- 
matics of gas within the bubble using various trac- 
ers. The velocity of particular clouds are reported 
as much as 20 km s _1 relative to the LSR; though 
it is not clear how these are related to the motion of 
the radio scattering plasma which may extend along 
much of the line of sight. 

Since these studies show a range of possible ve- 
locities in the local ISM, we considered a model 
in which the scattering medium is uniformly dis- 
tributed between the Earth and the pulsar, but is 
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Fig. 6. — Scintillation time scale for observation in 
April 1996 (solid lines) and in October 1998 (dashed 
lines), calculated for a uniform scattering medium 
with a uniform velocity (V a , Vs) relative to the LSR. 
The time scales have been reduced by the factor 1.20 
to account for the incipient strong scattering, as de- 
scribed in the text. The hashed region is where the 
computed time scales agree with the observations 
within ±cr in both years. 

moving with respect to the LSR at constant veloc- 
ity (V a ,Vg). We performed an exhaustive search in 
velocity space for a single medium velocity that si- 
multaneously matches the time scales 1996 and 1998 
data. Here the measured time scales were obtained 
from the thin screen models fitted to the structure 
functions, since these models provided a good tem- 
plate for the observations. 

We present the results in Figure ^| as contours 
of the predicted time scales (again reduced by the 
factor 1.2) in (V a ,Vs) for the two observing peri- 
ods. The hashed region shows where the calculated 
time scales lie within ±<r of the measurements. Its 
extreme boundaries are — 8 < V a < 1 km s _1 and 
— 7 < Vg < 4 km s _1 and just include the LSR 
itself. The worst case errors due to uncertainties 
in the pulsar proper motion are shown by the dot- 
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ted contours, which increases the velocity range to 
-12 < V a < 5 km s" 1 and -16 < V s < 7 km s _1 . 
Since these are worst case limits the associated con- 
fidence level is 90%. It would clearly be interest- 
ing to tighten the proper motion estimate. We em- 
phasize that this velocity only estimates the plasma 
motion in a plane normal to the direction toward 
the pulsar (at I = 140°, b = 31°) - plasma velocity 
along this direction is not constrained. Upon reflec- 
tion, it seems remarkable that the simplest uniform 
Kolmogorov model (with only one parameter - the 
integrated strength of turbulence) provides a satis- 
factory match to both observations. We also note 
that the (projected) solar velocity in this frame lies 
nearly 20 km s _1 away and is clearly excluded. 

To summarize, in the absence of other informa- 
tion about the nature of the scattering medium, we 
have two models that explain the observations. A 
uniform scattering medium which moves at nearly 
the LSR velocity and a local scattering region at a 
distance of about 195 pc which is stationary in the 
LSR. It is also evident that localized scattering re- 
gions at other distances would explain the observa- 
tions if they were moving at the right velocity. We 
examine the tradeoff between distance and veloc- 
ity in the context of other astronomical information 
about the local scattering environment in section 4. 

3.3. Spectra with an Inner Scale 

The possibility that the density spectrum is Kol- 
mogorov but includes an "inner scale" at which 
damping becomes important, has been considered 
by several investigators (e.g. Coles et al. 1987, Ro- 
mani et al. 1986, Spangler and Gwinn, 1990). The 
inner scale produces a rapid cutoff in the density 
spectrum which can be modelled with a gaussian 
multiplier of the form exp[— k 2 s 2 /4]. The inten- 
sity spectrum is directly related to the density spec- 
trum and shows the same multiplier. However, the 
structure function shows the effect much less dra- 
matically - it changes exponent from D(s) cx s 2 to 
D(s) cx s 5 / 3 above s = Sj if Sj <C r/. We also note 
that both the screen and extended medium curves 
in Figure|| approach the asymptotic 5/3 slope quite 
slowly; even the relatively steep thin screen model 
has only reached a slope of 1.44 at one tenth of 
the saturation level. Thus one should be quite cau- 
tious about interpreting the logarithmic slope of the 
structure functions at small lags as estimates of a. 

We used fits to the power spectra, rather than the 
structure function, to bound the inner scale. Wc 



tested both the thin screen and uniform medium 
models described earlier. The best fit theoretical 
spectra for thin screen and uniform medium models 
are shown in Figure || as solid and dashed lines re- 
spectively. In 1996 all models with Si < 6 x 10 6 m fit 
equally well, whereas in 1998 the fits at Si = 2 x 10 7 
m are significantly better. Although it is tempting 
to believe that an inner scale was detected in 1998, 
the detection depends on one or two points at the 
highest frequencies in the observed spectrum. These 
points are sensitive to the correction for system 
noise and intrinsic variation. We have assumed that 
both these corrections are white. A small amount of 
narrow band noise could significantly bias the high- 
est frequency spectral estimates, increasing or de- 
creasing them, depending on whether the noise is 
above or below the data range. Thus we have con- 
fidence only that the data indicate that Si < 2 x 10 7 
m. 

The only astrophysical plasma with which we can 
make a comparison is the solar wind. There an inner 
scale has been well determined (Coles and Harmon 
1989), and appears to be very close to the ion iner- 
tial scale, i.e. the Alfven speed divided by the ion 
cyclotron frequency. Thus Sj = 684 (N e cm -3 ) -1 / 2 
km depends only on the density. At a typical density 
of N e = 0.03 cm -3 we have s, = 4x 10 6 m. Spangler 
and Gwinn (1990) made a similar suggestion for the 
dissipation scale in the interstellar plasma, and pro- 
vided some evidence for inner scales <10 6 m for very 
heavily scattered lines of sight. It seems likely that 
the scattering takes place in regions of somewhat 
higher density than the mean interstellar medium, 
so it is likely that the inner scale is less than 4 x 10 6 
m. This is consistent with our observations but the 
bound is not very tight. The dynamic range of the 
observations would have to be increased by about 
an order of magnitude to provide a more interest- 
ing bound. An important consequence of the good 
spectral fits is to provide new and independent ev- 
idence in suport of the Kolmogorov power law for 
the spectrum of plasma density in the local ISM; 
this has become the canonical model, but is not of- 
ten tested critically. 

3.4. Scintillation Index 

In the foregoing analysis we have concentrated 
on interpretating the temporal variations charac- 
terized by the normalized structure functions. In 
weak scintillation the time scale depends only on 
the relative distribution of turbulence along the line 
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of sight, and is independent of its absolute level. 
We now estimate the absolute level of turbulence 
from the measured scintillation index and compare 
it with other observations of the same pulsar and 
also with models for the ISM electron density. The 
models describe the distribution of the variance of 
densi ty, pa rameterized by the constant in equa- 
tion ( A12 ) for the density spectrum. But, in gen- 
eral, the observations can only estimate the integral 
of along the line of sight, termed the scattering 
measure SM. 

Observations of ISS have shown extreme variabil- 
ity in SM between even rather closely aligned sight 
lines Previous workers (e.g. Cordes, Weisberg and 
Boriakoff, 1993 [CWB]; TC93) have modelled the 
distribution of as a locally uniform background, 
with a stronger and more clumpy component super- 
imposed. These models have been based primarily 
on pulsar measurements of the diffractive decorre- 
lation bandwidth 5vd in strong scattering. Though 
this parameter is easy to estimate over 1-2 hours, 
the repeated 327 MHz observations of Bhat et al. 
(1999a) show it to be surprisingly variable. For 
many of the 20 pulsars observed within 1 kpc, bvd 
varied by a factor of 5 or more over days to weeks. 
This variation is almost certainly caused by refrac- 
tive scintillation, which modulates the diffractive 
spectrum and is caused by much larger scales in the 
ISM than those responsible for the diffractive scintil- 
lation itself. As noted by Cordes, Pidwerbetsky and 
Lovelace (1986), and by Gupta et al. (1994), the re- 
fractive modulations also tend to bias the apparent 
5vd downwards, suggesting that the larger values of 
6i>d more accurately estimate the diffracting irregu- 
larities. With this preamble, we now compare SM 
derived from our scintillation index measurements 
with estimates based on 5ud measured for the same 
pulsar and also with SM from the TC93 and Bhat 
et al. (1998) models. 

Our observations of the scintillation index give 
SM ~ 3x 10~ 6 m~ 20 / 3 kpc, depending slightly on 
whether we assume a screen or uniform scattering 
medium (see Table 1). Also listed are SM esti- 
mated from two Svd observations (at 151 MHz by 
Rickett, 1970 and at 360 MHz by Cordes, 1986). We 
assume a uniform scattering medium in estimating 
SM from these measurements, and we noticed that 
Cordes (1986) used a formula higher by a factor 3.2 
than that in TC93. As discussed in Appendix C, we 
use the formula in TC93, since it agrees with the in- 
dependent analysis of Lambert and Rickett (1999). 



Even with this correction, the 151 MHz observation 
gives an SM 3.6 times higher than ours and the 
360 MHz one is 13 times higher. The two measure- 
ments of 5i>d differ by a factor 5 when scaled to a 
common frequency. From the discussion above, we 
may assume that the larger value of Svd, i.e. the 
smaller value of SM, better represents the under- 
lying density spectrum responsible for small scale 
ISS; but even then we have a factor of 3.5 discrep- 
ancy in SM. Both methods scale in the same way 
with screen distance and with pulsar distance, and 
so they cannot be reconciled by adjusting the dis- 
tances in either the screen or uniform scattering ge- 
ometries. In agreement with Bhat et al. (1999b), 
we conclude that we need a better understanding of 
the variability in 5vd observations; numerical sim- 
ulations will be necessary, since the existing the- 
ory assumes asymptotically strong scattering which 
is not strictly applicable on the shorter paths to 
nearby pulsars. 

We also calculated SM from the models of TC93 
and Bhat et al. (1998). The former is 23 times and 
the latter is 18 times our estimate. TC93 model 
the plasma as uniform turbulence in the Galactic 
disk of 700 pc scale height in Cjj. The recent 
model proposed by Bhat et al. (1998) adds detail 
near the Sun of a bubble of low turbulence sur- 
rounded by enhanced turbulence in a shell at 50-150 
pc. They propose a range of ellipsoidal models for 
the shell, which the line-of-sight to PSR B0809+74 
(I = 140°, b = 31.6°) crosses at 72 ± 13 pc from 
the Earth. The model specifies a low Cjy inside the 
bubble, a value of SM for the shell and the same 
C% as TC93 outside the shell. For our 310 pc path 
their shell contribution to SM dominates with 70% 
and the outer region provides 30% and the inner 
region less than 1%. Their model is closer to be- 
ing a "screen" than a uniform extended scattering 
medium. Since we find 18 times lower SM than in 
their model, we conclude that the enhanced turbu- 
lent shell is not continuous and must have a hole 
toward pulsar B0809+74. We note that of the 20 
pulsars observed by Bhat et al., the closest to PSR 
B0809+74 was PSR B0823+26 at the same Galac- 
tic latitude but 50 degrees away in longitude. Thus 
their observations are not in conflict with our low 
SM toward pulsar B0809+74. 

The large SM discrepancies with the scattering 
models could be resolved by substantially reducing 
the pulsar distance. If L were reduced by a fac- 
tor 5.5 it would bring our SM into line with that 
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(Jbservation 


Frequency 


Reference 


Model 


SAl 




(GHz) 






( m -20/3 kpc ) 


m w = 0.8 


0.93 


this work 


screen at 195 pc 


2.3 x 10~ e 


m w = 0.8 


0.93 


this work 


extended med. 


3.0 x 10~ 6 


5v d = 0.95 MHz a 


0.36 


Cordes(1986) 


extended med. 


4.1 x 10~ 5 * 


5v d = 0.1 MHz b 


0.15 


Rickett(1970) 


extended med. 


1.1 x 10~ 5 * 


Theory 




TC93 


extended med. 


7.1 x 10~ 5 


Theory 




Bhat et al. (1998) 


shell of bubble 


5.5 x 10~ 5 * 



Note. — * SM obtained using the formulae in Appendix C (and in TC93) 
"average of observations at 0.36 and 0.41 GHz scaled to .36 GHz 

° Value obtained from measurement of bandwidth that reduces scintillation index to 0.5, divided by 10 (see CWB) 



from TC93 and bubble models. If that were the 
case, at 56 pc this would make PSR B0809+74 the 
closest pulsar to the Earth and the mean electron 
density 0.1 cm~ 3 , 3 to 5 times higher than other 
estimates. For example, PSR B0950+08 has a mea- 
sured parallax distance of 130 pc and a mean elec- 
tron density 0.025 cm~ 3 . Though the distance to 
B0809+74 may be uncertain by as much as 50%, it 
cannot be reduced by the factor 5.5, and we con- 
clude that the turbulence level on this line of sight 
is much lower than described by either of the mod- 
els. From our estimate of SM the effective aver- 
age C 2 N toward B0809+74 - lO^m" 20 / 3 . This is 
comparable to that in the inner cavity of the bub- 
ble model. It also agrees with that toward PSR 
B0950+08 (I = 229°, b = 44°) from the measure- 
ments of Phillips and Clegg (1992), after applying 
the correction factor from Appendix C. The low C'jq 
for both of these pulsars suggests a quite large region 
of low turbulence in the hemisphere away from the 
Galactic center. The direction to our pulsar looks 
out from the local arm through an inter-arm region, 
which we suggest may have a lower plasma density 
than that described by the TC93 model. The low 
density might also be related to the local bubble, ex- 
cept that the X-ray observations of Snowden (1990) 
indicate that the bubble terminates between 50 and 
100 pc towards PSR B0809+74. 

4. Discussion 

We have found that the observations agree with 
the predictions for a Kolmogorov density spectrum 
in either a uniformly extended scattering medium 
moving at close to the LSR velocity or in screen at 
rest in the LSR, located in the range 170 - 220 pc 
from the Earth. We now consider whether these re- 
sults are consistent with other ideas about the local 



scattering environment. 

In the local bubble model of Bhat et al. (1998), 
enhanced scattering is expected at 72 ± 13 pc from 
the Earth. This is marked on Figure |] by a rectan- 
gular bar and is clearly inconsistent with the screen 
solution near 195 pc. If the scattering region were 
located at 72 pc and co-moving with the LSR, the 
predicted time scale would be about 30% longer 
than observed in 1996 and twice as long as observed 
in 1998. Consequently, we also examined screen 
models with a fixed position and a free velocity V scr 
relative to the LSR. 

We have measured weak scintillation time scales 
(r W 96 and T w gs) for two epochs. If a; = z Q /L denotes 
the fractional distance to the screen, then the results 
of Appendix A give r w = 0.98 v /ALx(l - x)/27r/|(l- 
x)V + xV p — V SCT \. With V appropriate for the 
two epochs there are two equations, from which two 
possible solutions for V scr can be found at a given x. 
In Figure [7] these two velocity solutions are plotted 
as partial ellipses (solid and dashed lines), anno- 
tated with sample distances. The LSR velocity of 
a screen necessary to match the two observed time 
scales is thus given by a vector from the origin to 
the appropriate point on either curve. The lower ve- 
locity solutions are shown with sample error ellipses 
corresponding to the uncertainties in the time scale 
estimates. There are also errors due to uncertainty 
in the pulsar proper motion, given moving the point 
P by ±10 km s _1 in either coordinate. These map 
in a non-linear fashion into an error in V scri which is 
~ 10 km s _1 for screens near the pulsar, decreasing 
for screens nearer the Earth. 

A screen at the proposed bubble distance of 72 pc 
needs a particular velocity of V a = —13.8 ± 3, Vs = 
9.4 ±4 km s _1 transverse to the pulsar line of sight. 
While this is a reasonable interstellar velocity, the 
problem remains that our estimate of SM is 18 
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Fig. 7. — Velocities relative to the LSR, projected 
on the plane normal to the line of sight toward pul- 
sar B0809+74. The pulsar velocity is marked by P. 
The Earth's velocity is shown for the observations 
in 1996 and 1998, lying on the dotted ellipse of the 
Earth's annual track, centered on the solar velocity. 
For each observation the line of sight cuts through 
the medium at a velocity that lies on a line between 
the Earth and pulsar velocities, "ext" marks 1 a 
region of allowed velocity for a uniform scattering 
medium that best matches the data. Screen models 
give pairs of solutions for the screen velocity, which 
follow partial ellipses. One solution is drawn dashed 
since its larger velocities make it less likely. At 5 
screen distances error contours are drawn for the 
allowed screen velocity. The plot does not include 
errors due to uncertainties in the proper motion. 
Also shown are the projected velocities of two rel- 
evant interstellar clouds in the solar neighborhood 
as discussed in the text. 

times smaller than in the Bhat et al. model. The 
solution from section 3.1 with a screen near 195 pc 
requires the velocity V a — — 13 ± 4, Vs = — 18 ± 10 
km s _1 indicated by the elongated error ellipse, 
rather than being at rest in the LSR. We have also 
marked as a parallelogram the (modest) velocity 



needed for scattering uniformly extended along the 
entire line of sight. The pulsar B0809+74 happens 
to lie at high ecliptic latitude and so the projected 
earth velocity (shown by the dotted ellipse in Fig- 
ure |7|) samples a large region of velocity space. It is 
likely that we can further constrain the scattering 
region by a series of weak ISS measurements at 2-3 
month intervals over the year. 

Figure [?] is also useful for displaying other veloc- 
ity measurements of the local ISM, two of which we 
have projected onto our line of sight and plotted 
as vectors. The arrow marked LIC is the projected 
LSR velocity of the local interstellar cloud observed 
to be entering the heliosphere from its backscattered 
Lyman alpha (see Lallement, 1998). This cannot 
be related to the bubble solution with a compara- 
ble velocity, since the bubble model has no scatter- 
ing nearer than about 70 pc. Genova et al. (1998) 
report measurements of Na I absorption lines over 
a wide range of galactic longitudes. They identify 
an interstellar cloud with velocity of 13.8 km s~ x 
(toward I — 225°, b — 5.4°) covering a large solid 
angle including PSR B0809+74. The arrow marked 
"Cloud P" is the projected velocity of this interstel- 
lar cloud, which is comparable to that of our sample 
screen solution at 195 pc. Genova suggests (private 
communication 1999) that it may be a large warm 
diffuse cloud that extends 70 to 120 pc. Of course 
the other question in such a comparison is that the 
optical absorbing cloud may not be physically the 
same as the ionized region of enhanced turbulence. 
We also note that the likely interstellar plasma ve- 
locities are relatively small. The Alfven speed in the 
warm ionized medium, which is assumed to be re- 
sponsible for the scattering, is thought to be about 
10 km s _1 , similar to that of cloud motions. 

Britton et al. (1998) have compared the angular 
scattering of pulsars with their temporal broaden- 
ing and concluded that the scattering is often close 
to the pulsar. They suggest that it may be due to 
a nebula associated with each pulsar. This scenario 
cannot explain our observations of B0809+74, pri- 
marily because it is thought to be relatively old, of 
order T = 10 8 yrs. If the scattering takes place at a 
distance R from the pulsar then the velocity of the 
scattering nebula V ne b with respect to the pulsar 
must be less than R/T. This is a very tight bound 
and such low velocities would produce much longer 
time scales than we observe. For example if R = 10 
pc, then V ne b <0.1 km s _1 . However from Figure 7 
we can see that a screen at a distance L — R = 301 pc 
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from the Earth must be moving at 12 to 14 km s -1 
with respect to the pulsar to explain our measured 
time scales. For scattering near the pulsar the time 
scale can be approximated by = 0.82T/ykR, 
which is 250 hrs at R = 10 pc. Evidently a model 
of scattering near the pulsar would only be tenable 
only if the age of B0809+74 were less than 0.6 x 10 6 
yrs. 

Our results can also be used to predict the inter- 
stellar scintillation for extragalactic sources B0716- 
+714 and B0836+710, which are within 5° of B0809- 
+74. Evidently, the predictions depend on whether 
or not we extrapolate the low level of scattering that 
we measure over 0.3 kpc out over the entire 1 kpc 
pathlength of the TC93 model. 

For example B0716+714 is known to show rapid 
optical and radio variations (Wagner et al. 1996). 
The optical variations cannot be due to ISS and 
so imply a very small optical size; if this same size 
applied at cm wavelengths then it must also scin- 
tillate as an effective point source in the interstel- 
lar medium. Extrapolating our low Cjy over the 
entire pathlength of the TC93 model, we predict 
m w ~ 0.4 , and spatial scale ~7x 10 8 m at 5 GHz. 
Because the transverse LSR velocity of the earth 
varies considerably in this direction, the associated 
time scale t w would vary from 50 hrs in October 
to 5 hrs in April. With the higher level in the 
TC93 model the point source m w ^2.1, indicating 
strong scintillation, which would split into diffrac- 
tive and refractive components; however, their spa- 
tial and temporal scales would not be very different 
from the estimates above. A more realistic source 
model would include a maximum brightness (mini- 
mum angular size), which would reduce the rms flux 
variation and increase the time scales; nevertheless 
the strong annual variability in time scale would re- 
main, if the scattering plasma does indeed share the 
motion of the LSR out to 1 kpc. 

The source B0836+710 is known to have a very 
compact component from 20 GHz VLBI observa- 
tions (Otterbein et al., 1998). If the diameter ob- 
served at 20 GHz applied at 5 GHz, it should also 
show interstellar scintillations, which are not ob- 
served (Quirrenbach, et al., 1992). However, at 5 
GHz the compact component appears to be self- 
absorbed, so predicting the expected scintillation 
index versus frequency will require a full model for 
the source structure versus frequency. 



5. Conclusions 

We have developed a method for observing mul- 
tiple epochs of weak scintillation of nearby pulsars 
whose proper motion has been measured. Our ob- 
servations for pulsar B0809+74 can be explained 
by a uniform distribution of plasma "turbulence" 
which is stationary in the LSR. The density spec- 
trum follows Kolmogorov law for scales between 
2 x 10 7 m and 10 9 m, i.e. the turbulence may have 
an inner scale <2 x 10 m. This bound is consis- 
tent with the hypothesis that the inner scale is the 
ion inertial scale, which would be about 4 x 10 6 m 
for an electron density of 0.03 cm~ 3 . The velocity 
of the scattering medium that best fits the obser- 
vations is in the range — 12 < V a < 5 km s _1 and 
— 16 < Vg < 7 km s -1 . The error bounds, which are 
dominated by error in the pulsar proper motion, in- 
clude the LSR. 

The scintillation is also consistent with scatter- 
ing from a localized region which is not stationary 
in the LSR, but this appears to be less likely for two 
reasons. First the scattering is considerably weaker 
than predicted by models of the distribution of tur- 
bulence in the local neighbourhood, so it is unlikely 
that the line of sight to B0809+74 is dominated by 
a local intensely scattering region. Second, if the 
scattering is dominated by a local region, it must 
be moving at a substantial velocity with respect to 
the LSR as shown in Figure^. While this is possible, 
it appears to be somewhat less likely. The scintilla- 
tion is not consistent with scattering in a "nebula" 
associated with the pulsar because the time scale 
of the scintillation is too short by several orders of 
magnitude. 

The authors would like to acknowledge the help 
of Michael Kramer who pointed out that pulsar 
B0809+74 has a substantial polarization. This lead 
to a reanalysis that significant changed the conclu- 
sions and improved this paper. We thank the direc- 
tor and Staff of EISCAT for the use of the Sodankyla 
Telescope and we thank the NSF for partial support 
under grant AST 94-14144. 

A. Weak Scintillation Theory 

We give here the theory of weak interstellar scin- 
tillation of a point source (pulsar) moving with re- 
spect to a scattering medium observed by a moving 
observer. Since the scattering is weak we can calcu- 
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late the contribution of a single thin screen and use 
the Born approximation to add the effect of many 
thin screens linearly. 

Consider a thin phase screen between a point 
source P and an observer O. The screen is at dis- 
tance Zp from P and z from O with z +z p = L. The 
screen introduces a phase change <j>(s{) at transverse 
coordinate s\. Then the Fresnel diffraction integral 
for the complex field at an observer (0,0, L) can be 
written: 

jk 



/(0,0,L) 



2ttz, 



exp[— jkL] J d 2 S\ 



exp 



.)<. 



(Al) 



where z e is defined by l/z e = l/z a + \jz v and k is 
the radio frequency propagation constant. The fac- 
tor jk/2nz e ensures that the field has unit average 
intensity (< //* >) at z = L. 

If now we consider the effect of moving the source 
to a transverse coordinate s p and the observer to s Q , 
which are small compared to z p and z Q , the field can 
be expressed as: 



f(s p , s Q , L) 
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(A2) 



Here 



$(s p ,s ) = k\s p - s a \ 2 /2L 
*off = s p z / L + s z p /L, 

where s c s is the transverse coordinate where a 
straight line from P to O intersects the screen. 



(A3) 



From eq (A2) we can find the intensity I(s p ,- 
s ,L) = ff* as a double integral. We need the 
intensity covariance 

Rj (tr p , <r , L) = ( I{s p , s Q , L) x 

I(s p + o-p, s + tr 0) L) ) (A4) 

Substituting the double integral form of I into the 
above equation generates the ensemble average of 
a quadruple integral over transverse screen vari- 
ables (si, «2, S3, S4). After transforming to variables 
which are the various sums and differences of these 
and performing two integrations we obtain the spa- 
tial covariance 
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Rj(o- p) er , L) 



exp 



-0.5Di(r,q) 



(2irz e y 
_ 3_k_ 

Ze 



d 2 r J d 2 q 
(r q + cr cff ■ q) 



(A5) 



where r = (si + S2 — S3 — S4)/2 and q = (sj — S2 + 
s 3 - s 4 )/2 and D 4 {r, q) = D(r) + D(r) - 0.5D(r + 
q) — 0.5D(r — q) and <r e ff is defined in an equation 
analogous to A3. Here the function D(r) charac- 



terizes the phase screen by the structure function of 
its phase at offset r. Standard analysis of a plane 
wave normally incident on a phase screen gives the 
correlation function for intensity at spatial offset s 
at a distance z from the screen (e.g. equation ( 5.3) 
of Tatarskii and Zavorotnyi, 1980). Equation ( A5) 
is related to this simpler result by replacing s by 
er c ff and z by z e . A general treatment for waves 
scattered by such a screen requires a separate anal- 
ysis for weak and strong scintillation. Under weak 
scintillation the important regions in r and q that 
dominate the integrations are where -D4 is small 
compared to one, and exp[— O.5D4] can then be ex- 
panded to first order and integrated over r. The 
result is often called the Born approximation [e.g 
equation (4.11) of Prokhorov et al. (1975)]. 

We are concerned with temporal variations caused 
as the observer moves at velocity V and the pul- 
sar moves at V p , both velocities with respect to the 
screen. The observer measures a temporal correla- 
tion function at time offset r, which is given by: 

R\ (r) =R s i (<t p = V p t, <t = V t, L) ( A6) 

The associated structure function is 

D\(t) = 2RU0) - 2RUr) (A7) 

The weak scintillation Born approximation gives: 



R\{r) 



d 2 K P^(k) 4sin 2 (K 2 z e /2fc) x 

exp[-j« ■ V cS t] (A8) 



Here P^k) is the spectrum of the phase introduced 
by the screen (ie for a normally incident plane wave) 
at wavenumber k, which corresponds to (kq/z e ) in 
the notation used above, and 



V cS = V Zp/L + Vpz /L 



(A9) 



which is the effective velocity of the point where a 
straight line from the pulsar to the observer inter- 
sects the screen. Note also that the low- wavenumber 
approximation in strong scattering is also given by 
this result with an extra exponential factor inside 
the integral exp[—D(nz e /k)]. For a layer of thick- 
ness A, we replace P<f,(n) by P^(k)A, where 

P'^k) = 2irr 2 e \ 2 P Nc {K X}Ky ,K z = 0,z o ) (A10) 
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where X — 2n/k and Pn c {k z Q ) is the three- 

dimensional power spectrum of the electron density 
in the layer at a distance z from the observer. If 
P^(k) is isotropic we obtain 

/•OO 

D){t) = 32ir 2 r 2 e \ 2 A / kcLk P Nc (k, z q ) x 
Jo 

sin 2 ( K 2 z e /2fc) [1 - J (Kaffir)] . (All) 

This equation is used to compute models for com- 
parison with the observations in section 3. With a 
simple power law density spectrum, 



Pjv.(k, z )=C 2 N {z )^ 
13=11/3, a = (3-2, 



(A12) 



the observed temporal structure function can be ex- 
pressed as 



D\(t, z ) = 2ml dscr(p) 



(A13) 



where 



d S cr(p) = 9a I 

Jo 



-1-a „;„2/,2 



sm z ( y y8) x 



P = 



9a = 



[1 - J (yp/2)] d y 
TV cS ^k/z e (z ) 
a2 1+a ' 2 T{l/2 + a/A) 
V^r(l - a/4) 



(A14) 



With these definitions d(oo) = 1, and for the Kol- 
mogorov spectrum (a = 5/3), we find (i scr (0.99) = 
0.5 and so the scintillation time scale is r w = 
0.98y/ Ze/k/V e s, and g a — 7.346. We also give ex- 
pressions for the scintillation index in terms of the 
screen coherence scale s scr as 



m 2 „ = 



ira/4 



(A15) 



r(l - a/2) sin(7ra/4) 
2y/z e (z )/k 

8Tr 2 r 2 e X 2 C 2 N (z )A x 
T(l - a/2) 
a2 a T(l + a/2) 



Here the coherence scale for the screen s scr is defined 
as the lateral spatial offset over which there is an 
rms difference of 1 radian in the screen phase. 

For an extended medium in weak scintillation we 
find -Dj(t) by simply integrating along the line of 
sight, where the screen thickness A becomes dz a . 



We proceed by assuming the medium to be "frozen" 
and refer V and V p to the rest frame for the 
medium. In completing the integral, we assume that 
Cjj is independent of z Q but include the z a depen- 
dence of the other variables: 

Dj(t, V p , r, 6 p ) = [ D\{t, z )dz /A 
Jo 

z e (z ) = z (L- z )/L (A16) 
V cS (z ) = V p \z 2 + 2z (L - z Q )r cos 9 p + 
{L-z ) 2 r 2 Y' 2 /L 

The net effect is that each screen of thickness dz 
contributes with a different time scale to the total 
signal. This rounds the abrupt transition to satura- 
tion that is found in the thin screen. The variation 
of V e f{ is determined by r = V a /V p and 8 p the angle 
between V p and V . Thus the detailed shape of the 
resulting structure function depends on the param- 
eters r and 9 p . At one extreme equal and parallel 
velocities (r = 1, P = 0), V c g is independent of dis- 
tance and the structure function is just a stretched 
version of d scr - In the case of anti-parallel velocities 
V e tY covers a wide range since there will be a distance 
at which it goes to zero and the corresponding time 
scale t w — > oo. In our application we know the pul- 
sar and observer velocities with respect to the Sun. 
Our basic model assumes that the medium (either 
screen or uniform extended medium) moves with the 
local standard of rest. But the theory can equally 
be applied to models with the medium in uniform 
motion, if V p and V represent the velocities with 
respect to the medium. 



B. Error in Structure Function Estimates 

We need an expression for the errors expected 
in estimates of the structure function, in terms of 
the time lag r and duration of the data sequence T 
and the characteristic time constant for the process 
t . We present here an analytical solution and a 
simulation. The estimator used can be written in 
continuous (symmetrical) form as 



d(t) = [ 1 Ax{t,T) 2 dt/r 

J-T'/2 



(Bl) 



where AX(t,r) = X(t - r/2) - X(t + r/2) and 
T" = T — t. This estimator is unbiased as: 



{D(r)) = {AX(t,r) 2 ) = D(r) 



(B2) 
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which is independent of t or T. The covariance in 
the estimators at time lags t\ and t 2 is then 



Cov[D( n ),D(T 2 )} = 



(B3) 



Substituting Bl and taking the ensemble average 
inside the integrals this yields 

Ti/2 f T^/2 



(B4) 



Cov[£>(ri),£>(r 2 )] = 

-T[/2 J-T'J1 
{AX(t 1 ,Tl) 2 AX(t 2 ,T 2 ) 2 ) 

dt x dt 2 /T[T^ - D{n)D{r 2 ), 



where T[ = T — \t\\ and a similar expression for T' 2 . 
We follow the method used by Jenkins and Watts 
(1969, pp 412-3, - JW) and assume that AX is a 
Gaussian random variable with zero mean. Using 
the well-known average of the product of four such 
variables, we obtain: 

Cov[D(n),D(T2)} = / / 

J-T[/2 J-T^/2 

2(AX(t 1 ,T 1 )AX(t 2 ,T 2 )) 2 dt 1 dt 2 /T[T/ ! . (B5) 

Following (JW), we change to sum and difference 
variables of integration. The region of integration 
becomes a parallelogram, over which the sum vari- 
able is directly integrable to give 



Cov[D(ti),£>(t 2 )] = 



T" 



4T(r, u, vf{T" - \r\)dr/(T" 2 - v 2 ) 



4F(r,u,v) 2 (l - \r\/v)dr/{T h 



(B6) 



where u = (|n| + |r 2 |)/2, v = \\r 2 \ - |n||/2 and 
T" = (T{ + T 2 ')/2 = T - u, and 

r(r, u, v) = 0.5[D(r - u) + D(r + u) - 

D(r - v) - D(r + v)]. (B7) 

Variance in structure function with power law be- 
havior 

The interesting special case of the variance in 
D(t) corresponds to v = u = r and T" — T' , for 
which the second integral in equation ( |B6| ) is zero. 
We proceed by considering the various asymptotic 
forms. For small lags, D increases (monotonically) 
with lag and saturates at = 2a 2 x . A character- 
istic time scale separates these regions and can be 



defined as D(r ) = 0.5-Dqo. There are limiting cases 
in which the observing sequence is long enough to 
include many independent variations (T » t„) and 
the opposite extreme (T«t„). 

First, consider T > r Q . Then with T > r > r D , 



the triangular weight in eq (B6) can be ignored and 
the upper limit of the integral can be extended to 
infinity: 

/>oo 

Var[£>(r)] - (4/T') / T(r,T,v = 0) 2 dr. (B8) 
Jo 

This gives the expected behavior for rms errors in 
D as proportional to \/t /T' . At large lags r > t q , 
D then approaches saturation and 



r(r, r, v = 0) ~ - D(r) = 2R(r) 



(B9) 



where R{r) is the covariance function is related 
to D(r) as indicated. When put back into equa- 
tion^ B8), the result is 4 times the variance in the 
conventional autocovariance estimator at zero lag 
i?(0), using the standard formula given by JW equa- 
tion (5.3.21) with lag equal zero. We can then see 
that the two methods for estimating a\ have the 
same rms error. Having discussed the errors for 
large lags, we now address small lags r <C t q '. 
Then T can be approximated as 



r(r, t,v = 0) ~ D"{r)r 2 /2. 



(BIO) 



Hence we see that for small lags the rms error in 
D(t) varies as r 2 , for any well-behaved D{t) func- 
tion. 

To make further progress, we assume a particular 
form for D{t). The cases of interest are where D is 
a power law cx r a for small lags, and we adopt the 
following simple model (with 1 < a < 2) 



D(T)=D 00 T a /(r: + T a ). 



(Bll) 



We can use this in equation(B8) for values of r > r, 
but the integral is not done simply. We have inves- 
tigated the integral numerically and find the follow- 
ing to be a useful approximation which connects the 
behavior for r < r G and r < r D . 



D^2 



(B12) 



The accuracy of this approximation is best judged 
from the simulations described below. The impor- 
tant point about this analysis is that for exponents 
a less than 2, as r decreases the error in D{t) de- 
creases more quickly than D(t) itself, and so the 
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fractional error in D decreases with small r, as il- 
lustrated in the simulations of Coles and Harmon 
(1989). 

The second limiting condition to be discussed 
is when the duration T is short compared to r . 



Time Series from one Simulation fo= 2.0 e*p=- 2.6667 



Then the triangular weight in equation(B6) must 
be included. For the variance in D we again set 
v = and since T <C t q we can approximate 
D(t) r~j Doo(T/T ) a . Approximating the integral 
for 2 > a > 1.5 gives 



DirU^D^T/T'fiT'/r.) 2 



(B13) 



This result when normalized by D(t) shows a frac- 
tional error varying cx (r/T') 2_a , again decreasing 
with small t. 

Eq. ( p36| ) also allows us to consider the covariance 
of D estimates at lags t% and t%. Examining the 



influence of non-zero v, equation(B7) shows that as 
v increases from zero it only has a substantial influ- 
ence when it is comparable to u. In other words the 
estimates are only independent when the difference 
in time lags is comparable to their mean. A conse- 
quence of this correlation between estimates of D is 
that one can estimate the exponent a with reason- 
able accuracy even with small time lags and data 
duration less than r Q . 

To investigate the validity of the foregoing anal- 
ysis and its approximations we have simulated ran- 
dom time series with a power law structure func- 
tion, applied the structure function estimator and 
studied its error properties. The process was spec- 
ified by its power spectrum which was P(f) — 
Po[{l + (//7o) 2 r (a+1)/2 - This has a high frequency 
form oc / — ° , and the corresponding structure 
function is cx r a for rf -C 1. 

We use the standard method of simulating the 
process in the frequency domain. Independent zero 
mean Gaussian random variables of variance P(/)/2 
are stored into the real and imaginary parts of the 
complex Fourier representation of the series. After 
making it an Hermitian spectrum, an FFT is used 
to obtain the (real) random variable X n at time 
samples t n . Figure B8(a) shows a single example 
with a = 5/3 and < n < 4095. The turnover 
frequency is set low, to give over two decades with 
power law behavior; this correspond to sample in- 
terval of r o /180 and T = 23r Q . We applied the 
structure function estimator eq (Q) to the time se- 
ries, and repeated the process on each new realiza- 
tion of the time series. This was repeated to ob- 
tain many independent estimates D{t), whose mean 
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Fig. B8. — Simulation of a time series with a power 
law structure function (oc i 5 / 3 ). (a) Single sam- 
ple of the time series (b) Ensemble average power 
spectrum, (c) Solid line is an estimate of the struc- 
ture function from the single realization shown in 
(a); dashed line is the theoretical average structure 
function and dotted lines are plus and minus the 
standard deviation in the structure function esti- 
mate from (d). (d) Cumulated estimate of the rms 
error in the structure function estimates from 100 
realizations (solid line) ; dashed line is the analytical 
approximation to the error given by equation (B12) 
dotted line shows the average structure function on 
the same scale, illustrating that the relative rms er- 
ror in the estimate improves at small time lags. 



and rms were computed versus r. In Figure B8(b), 
we show the ensemble average power spectrum; in 
Fig B8(c) we show as a heavy line an estimation D 
from the single time series of Fig B8(a); the dashed 
line is the theoretical average D(t); the dotted lines 
show this mean plus and minus two standard devia- 
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Normalized errors in D for fo= 2.0 exp=- 2.6667 




Fig. B9. — Simulation of a time series with a power 
law structure function (oc i 5 / 3 ). Set of 10 errors 
D(t) — D(t) normalized by the expected error from 
equation (B12). This shows that for r less than 
about 0.3 t the estimation errors are highly cor- 
related, and that for lags above t q the errors are 
correlated over a range in lag of about r . 

tions deri ved from 50 realizations of D. The fourth 
panel Fig |B8|(d) shows the cumulated estimates of 
the rms error in D(t) as a solid line and the ana- 
lytical approximation of equation( |B12] ) as a dashed 
line, which evidently provides a very useful approx- 
imation to the error. The dotted line shows D{t) 
with the logarithmic scales it illustrates how the 
fractional error in D becomes smaller for small r. 
The structure function estimates are plotted from 2 
to 2048 sample intervals. 

The degree of correlation between the D esti- 
mates at two time lags, can be found as a function 
of the two lags. Instead of displaying such a correla- 
tion surface, we simply plot in Figure B9 a set of 10 
errors D(t) — D{t) normalized by the expected er- 
ror from equation(B12). This shows that for r less 
than about 0.3 r G the estimation errors are highly 
correlated, and that for lags above r Q the errors are 
correlated over a range in lag of about r Q . 



C. Scattering Measure Estimates 

For a sc reen of thickness A at distance xL equa- 
tions (A16) together give the relation of to SM. 
With a = 5/3 this gives 



SM = C%A = 0.01465 x 
ml[x(l - x)L\/2ir]- 5 / 6 r- 2 \- 2 



For our measurements this gives SM ~ 2.8 x 
10 -6 m -20 / 3 kpc, for a screen at 234 pc and (by co- 
incidence) the same value for a screen at 72 pc (the 
edge of the bubble in the Bhat et al model). If, 
instead, we assume a uniform scattering medium, 
then 



SM = 0.0664 m 2 w (L\/2ir) 



-5/6„ 



(C2) 



and we obtain SM ~ 3.0 x 10 b m 2U ' 3 kpc for our 
pulsar measuremnt. 

SM estimates were also been found by Cordes 
(1986), in which he measured the diffractive decor- 
relation bandwidth Svj for various pulsars and con- 
verted to Cfj for an equivalent uniform scattering 
medium. However, the formula used, his equa- 
tion (6), derived by Cordes, Weisberg and Boriakoff 
(1993) [CWB] differs from those in the more recent 
model by TC93. For a spherical wave point source 
in a uniform scattering medium at distance L from 
the observer, the relation can be expressed in the 
general form based on CWB equation (6): 



SM = Av a+2 L-~ a ' 2 5i 



-a/2 



(C3) 



(CI) 



CWB derived A in their Appendix A and then in- 
creased it by a factor 6 to account for the effec- 
tive line-of-sight weighting. Cordes (1986) used the 
resulting expression in deriving an average Cfj = 
SM/L for ISS measurements of over 70 pulsars. 
However, in constructing their model for the distri- 
bution of dispersing and scattering electrons in the 
ISM, TC93 used a different formula (their equation 
8). It is expressed in terms of the scatter broaden- 
ing time, which they relate by t s — (2irSi>d) to 
the decorrelation bandwidth. Putting this into the 
same format as above, their formula corresponds to 
A = 6.3 x 10~ 4 where SM is in the units m _20 / 3 kpc, 
v is in GHz, L is in kpc and 8vd is in MHz. We note 
that this is about one third of the value (A = 0.002) 
used by CWB and Cordes (1986). In the study of 
scattering in the local bubble Bhat et al. (1998) used 
this larger value also. 

The work of Lambert and Rickett (1999) al- 
lows us to obtain an independent theoretical eval- 
uation for A. They computed the spherical wave 
two-frequency "diffractive" correlation function for 
three spectral models, including the Kolmogorov 
spectrum. Asymptotically strong scattering was 
assumed in a uniform scattering medium, leading 
to computed intensity correlation functions, from 
which precise values of the decorrelation bandwidth 
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were obtained. These are given in terms of a nor- 
malized decorrelation width v e i for each spectrum 
model; the value Vd = 2.320 was found for the Kolo- 
mogorov spectrum. The decorrelation bandwidth is 
then 5v<i — WdSQ S /(L\/2ir). Here sq s is the ob- 
served spherical wave field coherence scale, where 
the phase structure function equals one. It is related 
by standard formulae to SM (e.g. their equations 11 
and 13). Eliminating s 0s we find: 

_ (2^ rf )"/ 2 r(l + Q /2) a ( Q + l)2° 
cxt " 8^ 2 r(l - a/2)r e 2 C 2+«/2 ' ^ 

where c is the speed of light. With a = 5/3 and us- 
ing the same units as above, we find A = 6.25 x 10 -4 
in very close agreement with the formula used by 
TC93, and 0.315 times the value used in work fol- 
lowing Cordes (1986) and CWB. When comparing 
observations and theory with an accuracy of better 
than a factor of 3, observers should use the constant 



A = 6.3 x 10~ 4 in equation (C3). In particular, in 
the models of Bhat et al. should be reduced by a 
factor of three over their published values. We have 
made the appropriate changes in listing SM values 
for pulsar B0809+74 in Table 1. For a screen at 
distance xL from the Earth, a similar analysis leads 
to: 

A SC1 . = 0.179A cxt [x(l-a;)]- 5 / 6 (C5) 
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